Deep learning volumetric deformable registration

ABSTRACT

A method and system for automated deformable registration of an organ from medical images includes generating segmentations of the organ by processing a first and second series of images corresponding to different organ states using a first trained CNN. A second trained CNN processes the first and second series of images and the segmentations to deformably register the second series of images to the first series of images. The second trained CNN predicts a displacement field by minimizing a registration loss function, where the displacement field maximizes colocalization of the organ between the different states.

RELATED APPLICATIONS

This application claims the benefit of the priority of ProvisionalApplication No. 63/301,975, filed Jan. 21, 2022, which is incorporatedherein by reference in its entirety.

FIELD OF THE INVENTION

The present invention relates to a convolutional neural network(CNN)-based deformable registration algorithm to reduce computation timefor analysis of medical images such as CT and MRI.

BACKGROUND

Image registration is the process of identifying a spatialtransformation that maps two (pair-wise registration) or more(group-wise registration) images to a common coordinate frame so thatcorresponding anatomical structures are optimally aligned. In otherwords, a voxel-wise “correspondence” needs to be established between theimages. Automated co-registration is an essential feature of nearly allmedical imaging applications, however, co-registration can becomputationally expensive and particularly challenging for imaging ofdeformable organs, an extreme example of which is the lungs. Theregistration of inspirative and expirative lung CT images has importantmedical applications in the diagnosis and characterization of smallairway disease and emphysema. Using accurate pulmonary registration, thelocal lung ventilation can be reliably quantified. Lung registration isalso widely used in radiotherapy for estimation of tumor motion togetherwith localized ventilation measures, which can be included into planningto spare as much well-functioning tissue as possible from radiationdose.

Diseases affecting the small airways, such as chronic obstructivepulmonary disease (COPD),bronchiolitis obliterans related to stem celltransplantation or chronic lung allograft dysfunction, and cysticfibrosis, can manifest as pulmonary air trapping, which can goundetected on routine inspiratory chest CT. For each of these diseases,chronic inflammation and obstruction of the small airways limit the rateand volume of gas expulsion during expiratory phase, which is typicallydiagnosed via pulmonary function testing (PFT) through measurement ofFEV1 and FEV1/FVC. Although air trapping may sometimes be observed asmosaic attenuation on expiratory phase CT, evidence has shown thatdiffuse air trapping can be difficult to detect visually. In contrast,quantitative measurements on a dedicated inspiratory/expiratory lung CTprotocol can facilitate air trapping assessment, which has been shown toprognosticate both disease progression and mortality.

Several methods for quantifying air trapping have been proposed. Earlymethods approximated air trapping by measuring low attenuation areas(LAA) on expiratory phase CT images, however, these measurements can beconfounded by areas of emphysema. Other methods quantify air trapping byregistering inspiratory and expiratory phase images using lungdeformable registration, enabling regional discrimination of areas wheregas exchange is impaired by air trapping from areas that areemphysematous. These iterative deformable registration algorithmsincorporate diffeomorphic constraints to enforce transformationinvertibility (i.e., limit voxel “folding”) and improve registrationaccuracy, but require a significant amount of time, from minutes tohours, to complete, increasing computational cost and limitingfeasibility for use in routine clinical care.

Recently, deep convolutional neural networks (CNN) algorithms have shownpromise for performing deformable registration in a variety of organs,each with the potential to reduce computational time, while preservingaccuracy. With the advent of deep learning (DL), there have beensignificant advances in algorithmic performance for various computervision tasks, including medical image registration. Several CNNarchitectures have been proposed in recent years, each employingspecific architectural modifications to address the issue ofvanishing/exploding gradients that are common to deep networks, such asAlexNet, VGG, ResNet, and DenseNet. Among these, in medical imagesegmentation and registration, the most widely used architecture is theU-Net --an encoder-decoder style network with skip connections betweenthe encoding and decoding paths. The encoder contains severalconvolutional layers and pooling layers, which downsample the inputimage to a low resolution. The decoder includes deconvolution layerswith a matching number of layers to the encoder. Through the decoder,the feature maps are reconstructed to the original size of the inputimages. The U-Net utilizes several down- and up-sampling layers to learnfeatures at different resolutions, at the limited expense ofcomputational resources. It has been widely applied in various medicalimaging applications, including segmentation. The approach of thepresent invention seeks to exploit the benefits of such CNNarchitectures to address the challenges of accurately assessingpulmonary air trapping.

SUMMARY

A CNN-based deformable registration algorithm performs automatedco-registration of organs. Application of the algorithm to an input oflung CT image data provides accurate quantification of air trappingmeasurements in a broad cohort of patients with a wide range ofpulmonary air trapping. The inventive approach reduces inferenceruntime, improves lobar overlap, and reduces voxel “folding.” Thesefast, fully-automated CNN-based lung deformable registration algorithmscan facilitate translation of measurements into clinical practice,potentially improving the diagnosis and severity assessment of smallairway diseases.

In the disclosed embodiments, the inventive approach employs a CNN-basedalgorithm to perform deformable lung registration using CT image data.Accuracy of the approach is enhanced by incorporating mathematicalconstraints to ensure applicability for lobar air trappingquantification. Included among the constraints are one or more of: (a)lobar segmentations to preserve anatomic boundaries, (b) displacementfield regularization to encourage physically realistic transformations,and (c) the Jacobian determinant to limit nonanatomic voxel “folding.”For deformable registration of chest CT, the inventive CNN-basedalgorithm achieves greater lobar overlap, and faster runtime (418×) thanan iterative reference method.

In one aspect of the invention, a method for automated deformableregistration of an organ from medical images includes: receiving in acomputer processor configured for executing a trained convolutionalneural network (CNN) image data corresponding to a first series ofimages of the organ in a first state and a second series of images ofthe organ in at least one second state; propagating the first and secondseries of images using a first trained CNN to generate segmentations ofthe organ; propagating the first and second series of images and thesegmentations of the organ using a second trained CNN to deformablyregister the second series of images to the first series of images,wherein the second trained CNN is trained to predict a displacementfield by minimizing a registration loss function, wherein thedisplacement field is configured to maximize colocalization of the organbetween the first state and the at least one second state; andgenerating an output comprising a deformed image corresponding to the atleast one second state. The registration loss function may be acombination of one or more loss function selected from the groupconsisting of cross-correlation loss, displacement field loss, Diceloss, and Jacobian loss. In some embodiments, the image data comprisesvoxels and the first and second trained CNNs are 3D U-Net CNNs. Theorgan may be one or more lung, where the first state is inspiratory andthe second state is expiratory, and the image data comprises avolumetric lung CT scan. The method may further include extracting lungmeasurements from the segmentations and/or generating a disease map froma visual comparison of the first series of images and the deformedimage. The image data may be pre-processed by one or more of resizing toa standard resolution, scaling voxel attenuations, andaffine-registering the second series to the first series.

In another aspect of the invention, a system for automated deformableregistration of an organ from medical images includes at least onecomputer processor configured to:

acquire image data corresponding to a first series of images of theorgan in a first state and a second series of images of the organ in atleast one second state; execute a first convolutional neural network(CNN), wherein the first CNN trained to generate segmentations of theorgan using the first and second series of images; execute a second CNN,wherein the second CNN is trained to deformably register the secondseries of images to the first series of images using the first andsecond series of images and the segmentations of the organ, wherein thesecond CNN is further configured to predict a displacement field byminimizing a registration loss function, wherein the displacement fieldis configured to maximize colocalization of the organ between the firststate and the at least one second state; and generate an outputcomprising a deformed image corresponding to the at least one secondstate. The registration loss function may be a combination of one ormore loss function selected from the group consisting ofcross-correlation loss, displacement field loss, Dice loss, and Jacobianloss. In some embodiments, the image data comprises voxels and the firstand second trained CNNs are 3D U-Net CNNs. The organ may be one or morelung, where the first state is inspiratory and the second state isexpiratory, and the image data comprises a volumetric lung CT scan. Themethod may further include extracting lung measurements from thesegmentations and/or generating a disease map from a visual comparisonof the first series of images and the deformed image. The image data maybe pre-processed by one or more of resizing to a standard resolution,scaling voxel attenuations, and affine-registering the second series tothe first series.

BRIEF DESCRIPTION OF THE DRAWINGS

The patent or application file contains at least one drawing executed incolor. Copies of this patent or patent application publication withcolor drawing(s) will be provided by the Office upon request and paymentof the necessary fee.

FIG. 1 is a diagram showing the study design used to evaluate the lungdeformable registration algorithm (LungReg) using data from the publicCOPDGene dataset.

FIG. 2 is a flow diagram for the LungQuant system for fully automatedlung CT measurements.

FIG. 3A is diagram of the architecture for U-Net inspired LungSegalgorithm for lung and lobe segmentation; FIG. 3B is a diagram of thearchitecture for an implementation of the exemplary 3D U-Net CNN LungRegused to predict the displacement field defining the deformation; FIG. 3Cis a flow diagram of loss functions incorporated into the training ofthe lung deformable image registration algorithm, LungReg.

FIGS. 4A-4D plot Dice scores and cross-correlations across algorithmsand anatomic structures (FIGS. 4A, 4C) and paired differences with SyN(FIGS. 4B, 4D), respectively, for different lung regions.

FIG. 5 illustrates a case example comparing deformed images,segmentations, and displacement fields for LungReg_(α,β,γ) and SyN.

FIG. 6 shows ROC curves and AUCs for each lung deformable registrationalgorithm for predicting GOLD stages using %EM-LAA and %AT-ADM aspredictors.

DETAILED DESCRIPTION OF EMBODIMENTS

The examples and embodiments described herein are directed to lungregistration based on CT image data. It will be recognized by those inthe art that the inventive method and system are applicable toco-registration of other deformable organs and structures within a body,including, but not limited to, the circulatory, digestive, reproductive,and neurological systems. The method and system described herein aresimilarly not limited to registration of CT images but may also be usedin conj unction with other known imaging methods including, but notlimited to MRI and ultrasound.

An automated deep leaning-based approach for deformable registrationemploys a multichannel 3D convolutional neural network (CNN). Accordingto embodiments of the inventive approach disclosed herein, the CNN is aU-Net, a fully convolutional network that is known in the art, its namederived from its u-shaped architecture. (See, O. Ronneberger, et al.,“U-Net: Convolutional Networks for Biomedical Image Segmentation”,(2015) arXiv: 1505.04597 [cs.CV], incorporated herein by reference).Briefly, the U-Net employs repeated application of convolutions, eachfollowed by an activation function and a max pooling operation. Duringthe contraction, the spatial information is reduced while featureinformation is increased. The expansive pathway combines the feature andspatial information through a sequence of up-convolutions andconcatenations with high-resolution features from the contracting path.

The inventive CNN-based deformable registration approach was evaluatedusing the study design shown in FIG. 1 . A previously developed 3D lobarsegmentation CNN (LungSeg) was first applied to the 9,118 inspiratoryand expiratory series pairs from the COPDGene study, a multicenterobservational study designed to identify genetic factors associated withCOPD. The lung deformable registration algorithm (“LungReg”) was trainedto perform expiratory-to-inspiratory registration using the CT imagesand corresponding lobar segmentations. LungReg was then evaluated acrossseveral technical and clinical diagnostic metrics. Abbreviations used inthe figure and throughout the specification: “AUC” = area under thereceiver operating characteristic curve; “COPD” = chronic obstructivepulmonary disease; “DIR-Lab” = Deformable Image Registration Laboratory;and “GOLD” = Global Initiative for Chronic Obstructive Lung Disease.

A retrospective HIPAA-compliant study was approved by the institutionalreview boards of the participating institutions with waived requirementfor written informed consent. Study data included non-contrast CT andspirometric data for 9,649 patients collected between 2007-2012 obtainedas part of the COPDGene study. Exclusion criteria were missinginspiratory or expiratory CT (n = 531), resulting in 9,118 patients thatwere included in this study. Demographics are provided in Table 1 below.

TABLE 1 Characteristic Training Validation Testing All Demographics N7500 618 1000 9118 Men 53.27% (3995/7500) 51.78% (320/618) 52.70%(527/1000) 53.10% (4842/9118) Race (%NHW) 68.17% (5113/7500) 69.74%(431/618) 68.90% (689/1000) 68.36% (6233/9118) Age 59.65 ± 8.99 59.67 ±9.26 59.76 ± 9.00 59.67 ± 9.01 BMI 28.80 ± 6.22 29.38 ± 6.43 28.67 ±6.13 28.82 ± 6.23 Pack years 44.24 ± 24.91 44.73 ± 25.63 44.16 ± 23.9344.26 ± 24.85 Current smokers 51.81% (3886/7500) 51.29% (317/618) 52.00%(520/1000) 51.80% (4723/9118) Spirometry N 7453 615 993 9061 FEV1 2.26 ±0.91 2.30 ± 0.96 2.23 ± 0.94 2.25 ± 0.92 FEV1pp 76.83 ± 25.60 78.04 ±25.85 75.19 ± 25.30 76.73 ± 25.59 FVC 3.32 ± 1.00 3.34 ± 1.07 3.30 ±1.00 3.32 ± 1.00 FVCpp 87.27 ± 18.18 87.95 ± 19.27 86.44 ± 17.34 87.22 ±18.16 FEV1/FVC 0.67 ± 0.16 0.67 ± 0.16 0.66 ± 0.16 0.67 ± 0.16 GOLD 03286 284 413 3983 GOLD 1 581 46 74 701 GOLD 2 1419 109 203 1731 GOLD 3818 71 142 1031 GOLD 4 447 32 52 531 PRISm 902 73 109 1084 N 7500 6181000 9118 %EM-LAA 6.12 ± 9.78 5.79 ± 9.34 6.47 ± 9.87 6.13 ± 9.76 Perc15-915.85 ± 31.78 -915.07 ± 31.13 -916.14 ± 32.84 -915.83 ± 31.85 Meanatten. -838.44 ± 36.67 -836.90 ± 36.01 -838.68 ± 38.33 -838.36 ± 36.81Volume (liters) 5.47 ± 1.43 5.45 ± 1.41 5.47 ± 1.44 5.47 ± 1.43 N 7500618 1000 9118 %AT-LAA 22.72 ± 20.92 22.14 ± 19.46 23.94 ± 21.38 22.82 ±20.87 Perc15 -859 ± 59.23 -860.01 ± 56.34 -861.45 ± 60.09 -859.43 ±59.13 Mean atten. -735.31 ± 72.57 -734.47 ± 68.12 -738.23 ± 73.75-735.57 ± 72.41 Volume (liters) 3.23 ± 1.18 3.20 ± 1.09 3.26 ± 1.16 3.23± 1.17 Table notes: Categorical data are reported as a percentage andfrequency. Continuous data are reported as mean ± SD. Reported imagingstatistics are based on measurements produced by commercial software,Thirona, already included in the COPDGene dataset, used as a referencestandard for LungSeg performance evaluation. Characteristic definitions:%NHW =Percentage non-Hispanic White, BMI = Body mass index, FEV1 =Forced expiratory volume in 1 second, FEV1pp = Percent predicted FEV1,FVC = Forced vital capacity, FVCpp = Percent predicted FVC, GOLD =Global Initiative for Chronic Obstructive Lung Disease, PRISm =Preserved ratio impaired spirometry, %EM-LAA = Percent emphysema definedby low attenuation area on inspiratory, Perc15 = 15% percentile ofattenuation distribution, %AT-LAA = Percent air trapping defined by lowattenuation area on expiratory.

Inspiratory (200 mAs) and expiratory (50 mAs) CT series of the entirethorax were acquired without contrast using multidetector GeneralElectric, Philips, or Siemens scanners, each with at least 16 detectorchannels. Images were reconstructed using a standard soft tissue kernel,utilized submillimeter slice thickness (0.625-0.9 mm) and intervals(0.45-0.625 mm), with smooth and edge-enhancing algorithms.

The COPDGene dataset also included lung CT measurements computed bycommercialsoftware (Thirona): mean attenuation (Hounsfield units [HU]),lung volume (liters), attenuation 15th percentile (HU), and LAA-definedpercentage emphysema (% voxels ≤ -950 HU, %EM- LAA) or air trapping (%voxels ≤ -856 HU, %AT-LAA) from inspiratory and expiratory acquisitionsfor the lungs and each lung lobe. To assess localization error of theproposed lung deformable registration, we obtained 10inspiratory/expiratory CT series pairs publicly available through theDIR-Laboratory database. Each series pair contains 300 landmarkannotations uniformly distributed across the entirety of both lungs as areference standard.

Spirometric measurements included forced expiratory flow in one second(FEV1), percent predicted FEV1 (FEV1pp), and forced vital capacity (FVC)following administration of 180 µg of albuterol. Spirometricmeasurements were used to classify patients according to theGOLD stagingsystem for COPD severity.

Referring to FIG. 2 , to perform automated lung CT measurements, LungSegand LungReg can be combined together as “LungQuant”. Inspiratory andexpiratory series are first propagated through LungSeg, a CNN algorithmfor lung lobe segmentation. Expiratory series are then deformablyregistered to inspiratory series using the CNN-based lung deformableregistration algorithm LungReg to co-localize the lungs at the voxellevel. Lobar segmentations are then used to extract various lungmeasurements, including air trapping, from the inspiratory, expiratory,and deformed expiratory images. Inspiratory and deformed-expiratoryimages can also be used to create disease maps for visual assessment ofemphysema and air trapping severity. Details of the LungSeg and LungRegalgorithms are described below.

LungSeg: Referring to FIG. 3A, LungSeg is designed to segment the lungs,lung lobes, and trachea using a single volumetric lung CT as input. (SeeHasenstab et al. (“Automated CT Staging of Chronic Obstructive PulmonaryDisease Severity for Predicting Disease Progression and Mortality with aDeep Learning Convolutional Neural Network”, Radiol Cardiothorac Imaging2021;3(2):e200477, incorporated herein by reference)). LungSeg enablesquantification of lung CT measurements, including air trapping, acrossthe lungs and lung lobes.

The LungSeg architecture is designed to maximize the size and detail ofthe input array (192 × 192 × 192) while accounting for graphicsprocessing unit (GPU) memory limitations during training through rapiddown-sampling and up-sampling in the encoder and decoder layers of theU-Net. It may be noted that other architectures may work similarly wellfor the task of lobar segmentation. The goal of LungSeg is to enablequantification of lung CT measurements, including air trapping, acrossthe lungs and lung lobes.

LungSeg was developed using a retrospective convenience sample of 1037volumetric CT series from 888 current or former smokers undergoing lungcancer screening at multiple institutions (Inst. 1 & 2). Six hundredtwenty-one series were obtained as part of the National Lung ScreeningTrial and 416 series were sampled from another institution betweenAugust 2008 and October 2018. Demographic information was not availabledue to complete anonymization of imaging data. A summary of the scannersand imaging parameters (sample size, mean ± SD, range) used for LungSegdevelopment and testing are provided in Table 2 below.

TABLE 2 NLST Low Dose- Inst. 1 Low Dose - Inst. 1 Low Dose - Inst. 2 N416 341 280 Scanner Mfr. GE Medical 384 221 95 Philips 10 6 -- Siemens 782 185 Toshiba 15 32 -- Imaging Param. Tube Current (mA) 169.74 ± 227.24(40-1580) 92.68 ± 32.42 (50-225) 127.45 ± 60.95 (40-320) kVp 116.15 ±8.92 (100-140) 120.12 ± 1.53 (120-140) 120.93 ± 4.22 (120-140)Reconstruction Diameter (mm) 340.68 ± 38.70 (263-500) 342.34 ± 36.25(230-440) 345.21 ± 37.24 (260-468) Pixel Spacing (mm) 0.67 ± 0.08(0.51-0.98) 0.67 ± 0.07 (0.45-0.86) 0.62 ± 0.12 (0.34-0.91) SliceSpacing (mm) 1.27 ± 1.95 (0.3-10) 1.52 ± 0.36 (0.625-2) 1.85 ± 0.54(1-5) Slice Thickness (mm) 1.26 ± 1.85 (0.5-5) 2.07 ± 0.53 (1.0-3.2)2.08 ± 0.54 (1.25-5) Slice Dimensions (pixels) 512 × 512 512 × 512(433-768) × (512-768) Convolutional kernels I70f 2 - - B 3 - - B20f1 - - B30f - 80 184 B50f - - 1 B70f 4 2 - BONE 1 - - BONE+ 202 - - C 16 - CHST 2 - - E 1 - - FC03 2 - - FC08 1 - - FC10 - 32 - FC15 3 - - FC301 - - FC50 7 - - FC52 1 - - L 4 - - LUNG 20 - - SOFT 121 - - STANDARD 38221 95 YB 1 - -

Using Coreline, a commercially available software for artificialintelligence (AI)-based imaging diagnosis, ground truth segmentationswere created for each of the 1037 CT series. Segmentations includedmasks of the trachea and left lower, left upper, right lower, rightmiddle, and right upper lobes. Corrections to ground-truth segmentationswere made using an ITK-Snap, an open-source graphical application formedical image segmentation, to ensure accuracy

CT series and corresponding ground-truth masks were standardized tofeet-first supine orientation and resized to 192 × 192 × 192 resolutionusing cubic spline interpolation. CT voxel attenuations were then scaledby a factor of 1/3000 to stabilize training while maintaining theattenuation distribution across all series. Ground-truth masks werefirst represented by a single 192 × 192 × 192 array (S) with thefollowing encoding: 11 left lower, 12 left upper, 13 right lower, 14right middle, 15 right upper, 16 trachea. Edge gradients ||∇s(v)|| werethen calculated for each voxel v = (v_(x), v_(y), v_(z)) across theground-truth lobar masks to create masks of the fissures separating thelobes:

$\left\| {\nabla s(v)} \right\| = \left( {\frac{\partial s(v)}{\partial x},\frac{\partial s(v)}{\partial y},\frac{\partial s(v)}{\partial z}} \right),$

such that

$\frac{\partial s(v)}{\partial x} \approx s\left( {p_{x} + 1,p_{y},p_{z}} \right) - s\left( {p_{x},p_{y},p_{z}} \right).$

Similar approximations were used for

$\frac{\partial s(v)}{\partial y}$

and

$\frac{\partial s(v)}{\partial z}.$

Ground-truth masks were then converted to a 192 × 192 × 192 × 10 channelarray representing masks of the lungs, fissures, lung lobes, andtrachea. Referring to the right panel of FIG. 3A, convolutional blockscomprise a 3D convolutional, batch normalization (BN), and ReLU layer.

For characterization of the LungSeg CNN, image pre-processing was firstperformed. Inspiratory and expiratory series and segmentations werestandardized feet-first-supine and resized to 192 × 192 × 192 resolutionusing cubic spline interpolation. CT voxel attenuations were scaled by afactor of 1/3000. Expiratory whole lung masks were thenaffine-registered (translation, rotation, scaling, no shearing) toinspiratory whole lung masks for initial lung alignment using mutualinformation as an image similarity metric within the AdvancedNormalization Tools (ANTsPy v0.2.2) package, an optimized and validatedmedical imaging library for Python. Affine parameters were subsequentlyapplied to the expiratory series.

Following the flow shown in the left panel of FIG. 3A, LungSeg takes asinput a 192 × 192 × 192 × 1 array representing a volumetric lung CT andoutputs a 192 × 192 × 192 × 10 array representing the masks of thelungs, fissures, lung lobes, and trachea. Convolutional blocksthroughout the architecture comprise a 3D convolutional layer withkernel size 3, followed by a batch normalization layer with momentum0.9, and rectified linear unit (ReLU) layer. In addition to skipconnections, residual blocks are used throughout the architecture tomitigate vanishing gradients. The input image is initially processedusing a convolutional block with stride 2 and is followed by a stridedresidual block and max pooling layer, each further reducing spatialdimension by half; the strided residual block uses a stride of 2 in thefirst convolutional layer. Rapid downsampling in these layers allows thearchitecture to process subtle imaging features such as fissures thatare only visible at higher resolutions while reducing memoryconsumption. Feature maps are subsequently propagated through additionalresidual blocks and max pooling layers until arriving at 6 × 6 × 6 × 128resolution in the bottleneck. The decoder portion of the network thenuses a sequence of transposed convolutions and residual blocks toprocess and upsample the feature maps to input resolution for predictionof masks.

LungSeg was trained using the Dice loss function, equally weightedacross the 10 anatomic structures. Let S_(k) and ŝ_(k), where k = 1,..., 10 represent the ground-truth and predicted segmentations forstructure k. The Dice loss L_(seg) finds the average volumetric overlapbetween each of the respective structures,

$L_{\mspace{6mu} seg} = - \frac{1}{10}{\sum\limits_{k = 1}^{10}{2\frac{\left| {s_{k} \cap {\hat{s}}_{k}} \right|}{\left| s_{k} \right| + \left| {\hat{s}}_{k} \right|}}}.$

Dice values were equally weighted across structures to compensate forthe pixel imbalance between the sparse fissure masks and dense lungmasks, thus increasing focus on the fissure separations between thelobes.

CT series and ground-truth segmentations were partitioned at the patientlevel using a 90%/10% training-validation split (934 training, 103validation). LungSeg was then trained using the Adam stochasticoptimizer with an initial learning rate of 0.001 and a batch size of 2for 10 epochs. The pretrained model was then fine-tuned with a learningrate of 0.0001 while dynamically applying random rotations (-15 and 15degrees), shifts (-8 to 8 pixels), and zoom (80%-110%) to images duringtraining. Training was terminated when the validation loss stoppedimproving.

LungSeg was then applied to the 9,118 inspiratory and expiratory seriesfrom COPDGene for external testing. COPD CT measurements were calculatedfor each lung and lung lobe: 1) percentage lung voxels with lowattenuation area (LAA, attenuation ≤ -950) on inspiratory (%emphysema[%EM]), 2) percentage lung voxels with LAA (attenuation ≤ -856) onexpiratory (%air trapping LAA [%AT-LAA]), and for both inspiratory andexpiratory, 3) volume (liters), 4) attenuation 15th percentile (Perc15),and 5) mean attenuation (1). Agreement between LungSeg measurements andmeasurements from commercially available software included in COPDGenewas assessed.

Agreement for both inspiratory and expiratory measurements was verystrong (Intraclass Correlations (ICCs) >0.9) across each structure, withthe exception of mean attenuation for the right middle lobe. A lowermean attenuation ICC was also observed across all structures, likely dueto the inclusion of high attenuation vasculature within the lobarsegmentations used to derive the COPDGene measurements.

LungReg: For the LungReg CNN, imaging preprocessiong was performed usingthe same procedure used for LungSeg. The LungReg algorithm is based onthe VoxelMorph deformable registration framework, initially validatedfor brain registration. (See Balakrishnan G, et al., “VoxelMorph: ALearning Framework for Deformable Medical Image Registration”, IEEETrans Med Imaging 2019;38(8):1788-1800.) Inspiratory andaffine-registered expiratory images are first propagated through a 3DU-Net CNN, using the input-level fusion mode to spatial transformationfunction φ parameterized by displacement field u. The spatialtransformation is then applied to the affine-registered expiratoryimage, via spatial transformer layer, to deformably register theaffine-registered expiratory image to the inspiratory image. The U-Netis trained to predict a displacement field that maximizes colocalizationof anatomic structures between the inspiratory and deformed expiratoryimages resulting from the output level.

The CNN architecture of g_(w)(I,E), shown in FIG. 3B, is an expandedversion of the architecture reported by Balakrishnan, et al. The CNNtakes as input a 192 × 192 × 192 × 2 array representing I and Econcatenated along the channel axis. Output of g_(w)(I,E) is a 192 × 192× 192 × 3 array representing the predicted displacement field u. Theencoder of g_(w)(I,E) comprises sequences of 3D convolutions with stride2 and kernel size 3, each followed by a LeakyReLU layer with parameter0.2. The decoder alternates between convolutions, LeakyReLU layers, and3D upsampling. The spatial transformer E◦φ is a differentiable layerwithin the LungReg network that linearly interpolates across its eightneighboring voxel attenuations.

Gradient descent is used to optimize U-Net weights by minimizing a lossfunction comprising four components, which point to the U-Net since theyare used to optimize U-Net weights: (1) cross-correlation for imagesimilarity (Lcc), displacement regularization for smooth deformations(L_(φ)), (3) Dice overlap score for alignment of anatomic structures(L_(seg)), and (4) percentage of voxels with nonpositive Jacobiandeterminants (L_(jac)) to encourage transformation invertibility. Notethe segmentations are only used during LungReg training and are notrequired during inference time.

FIG. 3C provides a flow-diagram of loss functions incorporated into thetraining of LungReg, where back lines (“BK”) = forward propagation; bluelines (“BL”) = spatial transformations; orange lines (“OR”) = lossfunctions. Inspiratory (I) and affine-registered expiratory (E) imagesare propagated through a 3D U-Net CNN g_(w)(I,E) with weights w topredict a spatial transformation function φ parameterized bydisplacement field u. The spatial transformation is then applied to theexpiratory image (E◦φ) via spatial transformer layer to deform E suchthat attenuation values I(ν) and [I◦φ](ν) correspond to similar anatomiclocations. Given I, E, network g is trained to predict a displacementfield u that maximizes colocalization of anatomic structures between Iand E by minimizing a loss function, discussed below.

Cross-correlation loss (Lcc) encourages local similarity betweeninspiratory and affine-registered expiratory images while being robustto shifts in attenuation distribution attributed to higher density areasof the lungs typically observed in expiratory phase acquisitions. Let Brepresent a local volume of size b³ surrounding a voxel ν and let (Î)νrepresent inspiratory mean attenuation across voxels ν_(i) surrounding vsuch that

Î(v) = ∑_(v_(i) ∈ B)I(v_(i)); [Ê ∘ ϕ](v)

is defined similarly. The cross-correlation loss between I and E◦φ isdefined as

$\begin{array}{l}{L_{\mspace{6mu} CC}\left( {I,E,\phi} \right) =} \\{- {\sum\limits_{v \in \text{Ω}}\frac{\left( {\sum_{v_{i} \in B}{\left( {I\left( v_{i} \right) - \hat{I}(v)} \right)\left( {\left\lbrack {E \circ \phi} \right\rbrack\left( v_{i} \right) - \left\lbrack {\hat{E} \circ \phi} \right\rbrack(v)} \right)}} \right)^{2}}{\left( {\sum_{v_{i} \in B}\left( {I\left( v_{i} \right) - \hat{I}(v)} \right)^{2}} \right)\left( {\sum_{v_{i} \in B}\left( {\left\lbrack {E \circ \phi} \right\rbrack\left( v_{i} \right) - \left\lbrack {E \circ \phi} \right\rbrack(v)} \right)^{2}} \right) + \text{ε}}}.}\end{array}$

ε is a parameter designed to stabilize network training. Smaller valuesof Lcc imply higher cross-correlations between images, indicative ofstronger image similarity. Hyperparameters b and ε were selected duringtraining.

Displacement field loss (Lφ: Physically realistic transformations areencouraged by regularizing the smoothness of the displacement fieldthrough a constraint on the magnitude of spatial gradients within thedisplacement field using the following loss:

$L_{\mspace{6mu}\phi}(\phi) = {\sum\limits_{v \in \Omega}\left\| {\nabla u(v)^{2}} \right\|}.$

where

$\nabla u(v) = \left( {\frac{\partial u(v)}{\partial x},\frac{\partial u(v)}{\partial y},\frac{\partial u(v)}{\partial z}} \right)$

and

$\frac{\partial u(v)}{\partial x} \approx u\left( {v_{x} + 1,v_{y},v_{z}} \right) - u\left( {v_{x},v_{y},v_{z}} \right).$

Approximations for

$\frac{\partial u(v)}{\partial y}$

and

$\frac{\partial u(v)}{\partial z}$

are performed similarly.

Dice loss (Lseg): During LungReg training, overlap of the trachea andfive lung lobes corresponding to the inspiratory and deformed expiratorysegmentations is encouraged using a Dice loss function. The Dice loss iscalculated as the average Dice across the six structural segmentations:

$L_{\mspace{6mu} seg}\left( {s^{(I)},s^{(E)},\phi} \right) = - \frac{1}{6}{\sum\limits_{k = 1}^{6}{2\frac{\left| {s_{k}^{(I)} \cap \left( {s_{k}^{(E)} \circ \phi} \right)} \right|}{\left| s_{k}^{(I)} \right| + \left| {s_{k}^{(E)} \circ \phi} \right|}}}.$

Jacobian loss (Ljac): An additional displacement field penalty isimposed by regularizing its Jacobian determinants to encouragedeformations that result in fewer regions of noninvertibility (i.e.,“foldings”). Let J(ν) = ∇ϕ(ν) represent the Jacobian matrix of φ forvoxel v and let det[J_(φ)(ν)] represent its determinant. The Jacobianloss imposes a stronger loss penalty on deformation fields with a largenumber of nonpositive Jacobian determinants,

$L_{\mspace{6mu} jac} = \frac{1}{2}{\sum\limits_{v}\left( {\left| {\det\left\lbrack {J_{\phi}(v)} \right\rbrack} \right| - \det\left\lbrack {J_{\phi}(v)} \right\rbrack} \right)}$

The LungReg loss function is a linear combination of the four losses,

L_(reg) = L_(CC) + αL_(φ) + βL_(seg) + γL_(jac),

where α, β, γ > 0 are parameters that control the impact of therespective loss component on the LungReg loss function.

For training the LungReg CNN, series pairs and segmentations wererandomly partitioned to 7500 for training, 618 pairs for validation and1000 pairs for testing. LungReg was then trained using the Adamstochastic optimizer with an initial learning rate of 0.0001 and a batchsize of 1 for 20 epochs (150,000 steps). For each run, training wasterminated when the validation loss stopped improving. Since a singletraining run lasted ~48 hours, a greedy search across hyperparameterswas not feasible. The cross-correlation loss was assigned a window ofsize b = 9 and ε was set as 10 since smaller values of ε tended todestabilize training, possibly due to lack of attenuation variabilitywithin many of the cross-correlation windows. We optimized α, β, and γin successive order, first optimizing α: α = {0, 0.01, 0.1, 0.5, 1, 2};β = 0; γ = 0, then β: α = 1;β = {0,0.01,0.1,0.5,1,2};γ = 0, then γ:α =1;β = 0.1;γ = {0,10⁻⁶, 2 × 10⁻⁶, 10⁻³,10⁻¹,1}. We denote the final modelwith optimal regularization parameters, α= 1;β = 0.1;y = 2 × 10⁻⁶ as“LungReg_(α,β,γ).” Note that training of LungReg is completelyunsupervised when β = 0 and does not require a ground-truth. LungRegwhen β ≠ 0 is a semisupervised algorithm that requires the ground-truthlobar segmentations. However, ground-truth lobar segmentations are notrequired during test time. The CNN was trained in Python v3.6 using thetensorflow-gpu v2.2.0 deep learning library on a NVIDIA Quadro RTX 8000graphics card.

For testing, LungReg_(α,β,γ) was compared with affine registration as anaive benchmark, an iterative symmetric diffeomorphic registration (SyN)deformable registration algorithm and versions of LungReg withalternative loss functions excluding L_(seg) and L_(jac) (LungReg_(α))or L_(jac) (LungReg_(α,β)). Using the testing set, algorithms werecompared using Dice overlap between inspiratory and registeredexpiratory segmentations, cross-correlation between inspiratory andregistered expiratory series, and percentage voxels with nonpositiveJacobian determinant. Spatial accuracy for each algorithm was alsoassessed using landmark colocalization error (LCE) measured by 3DEuclidean distance in millimeters (mm) between inspiratory and deformedexpiratory landmarks for the 10 DIR-Laboratory reference cases.

CPU-based runtimes for affine, SyN, and LungReg inference on the testingset were recorded in seconds. GPU-based runtimes for LungReg inferencewere also recorded.

Performance Analysis

We computed attenuation difference air trapping (%AT-ADM), defined asthe percentage nonemphysematous voxels with attenuation differences ≤100 HU between the deformed expiratory and inspiratory series. Agreementbetween %AT-ADM measurements computedusing LungReg or SyN was assessed.We also compared the ability of LungReg and SyN to predict thespirometrically-defined GOLD stages using %EM-LAA and %AT-ADM aspredictors.

Statistical analysis was performed using RStudio (v3.6.1). Agreement wasassessed using intraclass correlation (ICC). Dice overlapscores,cross-correlations, and percentage nonpositive Jacobian determinantswere compared across algorithms using paired t tests. LCE for eachalgorithm was compared using a linear mixed effects model to account forcorrelations between within-patient observations, with landmarks nestedwithin patients as random effects and a four-level categorical variablerepresenting SyN and each LungReg algorithm as a fixed effect.Bonferroni correction was used to control for familywise error rate. 95%confidence intervals were analytically calculated as appropriate.Statistical significance was assessed using a 5% Type I error threshold;any use of the word “significance” refers to statistical significance.Runtimes were reported descriptively. GOLD stage prediction wasperformed using logistic regressions and assessed using receiveroperating characteristic (ROC) curve analysis; areas under the curve(AUCs) were calculated andcompared using bootstrapping.

Dice scores and cross-correlations for each algorithm and correspondingpaired differences with SyN are shown in FIGS. 4A-4B and 4C-4D,respectively. The x-axis labels in each plot are: Whole = both lungs; LL= left lung; RL = right lung; LLL = left lower lobe; LUL = left upperlobe; RLL= right lower lobe; RML = right middle lobe; RUL = right upperlobe; and TRA = trachea. α, β, γ correspond to LungReg with thecross-correlation loss, segmentation loss, and Jacobian loss,respectively. Significance and direction of paired differences areindicated by (+) and (-).

As expected, affine registration consistently produces significantlylower Dice scores and cross-correlations (P values < 0.001) across allalgorithms and structures. LungReg algorithms significantly outperformedSyN for all structures, except for the right middle and right upperlobes, where LungReg without the use of segmentations during training,and hence without any emphasis on lobar boundaries, producedsignificantly lower Dice (P values < 0.001). LungReg_(α,β) andLungReg_(α,β,γ) with training segmentations showed a significantincrease in Dice overlap (P values < 0.001). LungReg cross-correlationwas significantly greater than affine and SyN (P values < 0.001), butdifferences withSyN were unsubstantial in magnitude for the individuallung structures.

Percentage voxels with nonpositive Jacobian determinants for eachalgorithm is shown in Table 3, which lists the performance metricscomparing affine, SyN, and LungReg registration algorithms.

TABLE 3 Method Dice % Folding Voxels GPU Runtime (s) CPU Runtime (s)Affine 0.81 ± 0.04 \ \ 1.81 ± 0.21 SyN 0.93 ± 0.02 0.10 ± 0.13 \ 418.46± 246.47 LungReg_(α) 0.93 ± 0.03 0.51 ± 0.44 1.03 ± 0.04 2.27 ± 0.15LungReg_(α,β) 0.95 ± 0.02 0.50 ± 0.43 1.00 ± 0.03 2.28 ± 0.15LungReg_(α,β,γ) 0.95 ± 0.02 0.04 ± 0.05 0.99 ± 0.03 2.27 ± 0.15

Note that Dice scores are averaged across the five lung lobes. LungRegalgorithms incorporating segmentations into training showed improvementsin overlap of lung structures. We observed large percentages of“folding” voxels for LungReg_(α,β) and LungReg_(α,β,γ) relative to theother algorithms. However, incorporation of the Jacobian loss(LungReg_(α,β,γ)) reduced the percentage of folding voxels below SyNwhile maintaining lobar overlap. LungReg CPU-and GPU-runtimes duringinference are much faster than SyN.

LungReg algorithms without the Jacobian loss showed five times thepercentage of nonpositive Jacobian determinants than SyN. However,LungReg_(α,β,γ) with the Jacobian loss had a significantly lowerpercentage of these voxels (P value < 0.001) than all other algorithms,including SyN, and with less variability.

On CPU, affine runtime (see Table 3) was significantly faster than thedeformable registration algorithms (P values < 0.001). LungReg runtimeswere magnitudes faster (418×) than SyN (Pvalues < 0.001), withconsiderably less runtime variability. LungReg implementation onGPUreduced runtime by half, requiring only ~ 1 second on average forinference.

The landmark colocalization error (LCE) was determined for eachdeformable registration algorithm across the 10 DIR-Laboratory referencecases. The error (in mm) is provided in Table 4. Average LCE for eachalgorithm was 7.21 to 7.81 mm for LungReg and 6.93 mm for SyN. SyN hadthe lowest LCE across all algorithms (P value < 0.001), butLungReg_(α,β,γ) with the Jacobian loss had the lowest LCE of theCNN-based algorithms (P value < 0.001).

TABLE 4 SyN LungReg_(α) LungReg_(α,β) LungReg_(α),_(β,γ) Case 1 8.79 ±6.14 10.44 ± 5.86 10.71 ± 5.73 9.98 ± 5.39 Case 2 7.35 ± 5.24 9.92 ±5.95 10.30 ± 6.32 9.20 ± 5.52 Case 3 5.61 ± 2.89 4.24 ± 2.38 4.57 ± 2.774.21 ± 2.46 Case 4 7.47 ± 5.12 7.76 ± 4.39 7.89 ± 4.59 7.51 ± 4.40 Case5 4.88 ± 3.15 6.83 ± 3.44 7.49 ± 3.74 6.62 ± 3.38 Case 6 6.52 ± 4.416.66 ± 4.04 6.60 ± 3.98 6.09 ± 3.53 Case 7 5.49 ± 2.75 5.03 ± 2.95 5.13± 2.93 4.76 ± 2.77 Case 8 6.38 ± 3.84 7.15 ± 4.04 6.88 ± 3.99 6.82 ±3.96 Case 9 5.25 ± 3.23 6.31 ± 3.98 7.00 ± 4.87 6.54 ± 4.44 Case 1011.59 ± 5.94 10.02 ± 5.66 11.49 ± 6.13 10.39 ± 5.92 Avg 6.93 ± 4.27 7.44± 4.27 7.81 ± 4.50 7.21 ± 4.18

SyN LCEs were significantly lower than LungReg LCEs on average by 0.28to 0.88 mm (P value < 0.001). LungReg_(α,β,γ) with the segmentation losssignificantly increased LCE, relative to LungReg (P value < 0.001).However, incorporation of the Jacobian loss significantly reduced LCEbelow the other LungReg algorithms (P value < 0.001), suggesting thereduction of voxel “folding” improved registration accuracy.

Example 1

A case example comparing deformed images, segmentations, anddisplacement filed for LungReg_(α,β,γ) and SyN in a test set participantis provided in FIG. 5 . In the figure, dashed white lines are overlaidinspiratory segmentations. Positive displacements (red) areposterior-to-anterior, left-to-right, and inferior-to-superior. Deformedexpiratory images appear similar across LungReg_(α,β,γ) and SyN.However, overlap of lung structures improves, especially for the rightmiddle lobe (RML). Displacement fields suggest similar anatomictransformations between algorithms, but with greater emphasis on thelung boundaries, as evidenced by the lung outline presence visible ineach LungReg_(α,β,γ) field. Incorporation of masks into LungReg_(α,β,γ)training improves overlap of lung structures, especially for the rightmiddle lobe (LungReg Dice 0.93 versus SyN 0.87). Visual assessment ofdisplacement fields suggests similar anatomic transformations betweenalgorithms, but with greater emphasis on the lung boundaries.

Intraclass correlation coefficients (ICCs) between each LungRegalgorithm and SyN are listed in Table 5 below, indicating strongagreement (ICCs: 0.98-0.99) between LungReg and SyN air trappingmeasurements for each lung structure. Air trapping measurements appearrobust to the inclusion of the segmentation and Jacobian loss functions.FIG. 6 provides receiver operating characteristic (ROC) curves and areaunder the receiver operating characteristic curves (AUCs) for each lungdeformable registration algorithm for predicting Global Initiative forChronic Obstructive Lung Disease (GOLD) stages with use of percentageemphysema low attenuation area, or %EM-LAA, and percentage airtrapping-attenuation difference map, or %AT-ADM, as predictors.Algorithms showed near-identical performance for each respective GOLDstage (1-4). Areas under the receiver operating characteristic curvesfor the algorithm for performing deformable registration of lung CT(LungReg) were not significantly different from those for symmetricdiffeomorphic registration (iterative) (SyN), which suggests thatLungReg air trapping measurements could supplant SyN air trappingmeasurements.

TABLE 5 Structure LungReg_(α) LungReg_(α,β) LungReg_(α,β,γ) Whole 0.99[0.99, 0.99] 0.98 [0.98, 0.99] 0.98 [0.98, 0.99] LLL 0.98 [0.98, 0.99]0.98 [0.98, 0.98] 0.98 [0.98, 0.98] LUL 0.98 [0.98, 0.98] 0.98 [0.98,0.99] 0.98 [0.98, 0.98] RLL 0.99 [0.99, 0.99] 0.98 [0.98, 0.99] 0.99[0.98, 0.99] RML 0.99 [0.98, 0.99] 0.98 [0.98, 0.98] 0.98 [0.98, 0.99]RUL 0.99 [0.99, 0.99] 0.99 [0.98, 0.99] 0.99 [0.98, 0.99]

The inventive U-Net inspired CNN disclosed herein performs deformableregistration across the inspiratory and expiratory phases of therespiratory cycle to enable quantification of lobar air trapping. Ahybrid loss function that incorporates lobar segmentation overlap andthe Jacobian determinant improves algorithm performance, achieving animportant goal of reducing inference runtime from as much as ~15 minutesto ~2.25 seconds on CPU and ~ 1 second on GPU, without loss of accuracy.The final model, LungReg_(α,β,γ), achieved greater overlap of lungstructures and comparable spatial accuracy to the iterative SyNalgorithm, while generating displacement fields with fewer regions ofnonanatomic non-invertibility (“folding” voxels). Further, when weapplied LungReg_(α,β,γ) to CTs from phase 1 of the COPDGene study, thealgorithm exhibited similar ability to predict spirometric GOLD stage,which requires voxel-wise registration of inspiratory and expiratoryphase images.

The accuracy of the inventive algorithm is partially attributed to theincorporation of the lobar segmentation loss in the hybrid lossfunction. The addition of a segmentation loss encourages lobar overlapduring training, thereby improving colocalization of lung structures atthe lobar boundaries. As a result, we observed improved lobar overlap,particularly for the right middle lobe. This is particularly importantfor computing regional lung characteristics, such as air trapping, atthe lobar level.

Since the lung lobes span large territories, deformations within thelobes can be relatively unconstrained and subject to areas of voxel“folding” (i.e., areas of physically unrealistic and noninvertibletransformations measured by the percentage of nonpositive Jacobiandeterminants). To enforce transformation invertibility (i.e., to reducevoxel “folding”), some CNN deformable registration algorithms explicitlyrestrict transformations to be diffeomorphic. In contrast, the inventiveapproach uses a strategy that incorporates the Jacobian as an additionalloss component, providing improved landmark colocalization error whilemaintaining lobar overlap and reduced “folding” voxels to 0.04%,outperforming the diffeomorphic SyN algorithm. Overall, LungReg’sperformance, measured by Dice overlap (0.95), voxel “folding” (0.04%),landmark colocalization error (7.21 mm), and runtime (~ 1 second) wascomparable to other methods reported in the literature, reporting92%-93% Dice, 0.02 to 2.1% for voxel “folding”, 1.39 to 7.98 landmarkcolocalization error, and 0.2-2 seconds runtime. The inventive approachexpands on prior efforts, recognizing the need to address multipleendpoints for deformable registration, and ultimately confirms theutility of CNN-based deformable registration to quantify air trappingand stratify patients from the COPDGene study.

It may be noted that due to long training times, LungReg hyperparameterswere optimized successively (as described above). Alternative CNNarchitectures were not thoroughly explored for segmentation andregistration, however, it will be apparent to those of skill in that artthat other CNN architectures may be used. Efficient joint hyperparameteroptimization may provide further benefit as may additional/alternativeCNN architectures inthe context of lung measurements. In addition,segmentation overlap following deformable registration was calculatedusing lobar segmentations inferred by a previously developed lungsegmentation CNN, rather than manually hand-drawn 3D segmentations.Manual 3D segmentation of 1,000 series pairs would not be feasible suchthat the tradeoff of the high volume of test data outweighed thepotential benefits of a ground truth defined by manual segmentation.Since the test data included data solely from the COPDGene study, whichused a standardized imaging protocol, it may not reflect the variabilityof scans collected using clinical imaging acquisition protocols.

The inventive CNN-based deformable lung registration algorithm disclosedherein accurately quantifies air trapping measurements in a broad cohortof patients with a wide range of pulmonary air trapping, reducesinference runtime, improves lobar overlap, and reduces voxel“folding.”Fast, fully-automated CNN-based lung deformable registration algorithmscan facilitate translation of these measurements into clinical practice,potentially improving the diagnosis and severity assessment of smallairway diseases.

1. A method for automated deformable registration of an organ frommedical images, the method comprising: receiving in a computer processorconfigured for executing a trained convolutional neural network (CNN)image data corresponding to a first series of images of the organ in afirst state and a second series of images of the organ in at least onesecond state; propagating the first and second series of images of theorgan using a first trained CNN to generate segmentations of the organ;propagating the first and second series of images and the segmentationsof the organ using a second trained CNN to deformably register thesecond series of images to the first series of images, wherein thesecond trained CNN is trained to predict a displacement field byminimizing a registration loss function, wherein the displacement fieldis configured to maximize colocalization of the organ between the firststate and the at least one second state; and generating an outputcomprising a deformed image corresponding to the at least one secondstate.
 2. The method of claim 1, wherein the registration loss functionis a combination of one or more loss function selected from the groupconsisting of cross-correlation loss, displacement field loss, Diceloss, and Jacobian loss.
 3. The method of claim 1, wherein the imagedata comprises voxels.
 4. The method of claim 3, wherein the first andsecond trained CNNs are 3D U-Net CNNs.
 5. The method of claim 1, whereinthe organ comprises one or more lung, the first state is inspiratory andthe second state is expiratory, and the image data comprises avolumetric lung CT scan.
 6. The method of claim 5, further comprisingextracting lung measurements from the segmentations.
 7. The method ofclaim 5, further comprising generating a disease map from a visualcomparison of the first series of images and the deformed image.
 8. Themethod of claim 1, further comprising pre-processing the image data byone or more of resizing to a standard resolution, scaling voxelattenuations, and affine-registering the second series to the firstseries.
 9. A system for automated deformable registration of an organfrom medical images, comprising: at least one computer processorconfigured to: acquire image data corresponding to a first series ofimages of the organ in a first state and a second series of images ofthe organ in at least one second state; execute a first convolutionalneural network (CNN), wherein the first CNN trained to generatesegmentations of the organ using the first and second series of images;execute a second CNN, wherein the second CNN is trained to deformablyregister the second series of images to the first series of images usingthe first and second series of images and the segmentations of theorgan, wherein the second CNN is further configured to predict adisplacement field by minimizing a registration loss function, whereinthe displacement field is configured to maximize colocalization of theorgan between the first state and the at least one second state; andgenerate an output comprising a deformed image corresponding to the atleast one second state.
 10. The system of claim 9, wherein theregistration loss function is a combination of one or more loss functionselected from the group consisting of cross-correlation loss,displacement field loss, Dice loss, and Jacobian loss.
 11. The system ofclaim 10, wherein the image data comprises voxels.
 12. The system ofclaim 10, wherein the first and second CNNs are 3D U-Net CNNs.
 13. Thesystem of claim 9, wherein the organ comprises one or more lung, thefirst state is inspiratory and the second state is expiratory, andwherein the image data comprises a volumetric lung CT scan.
 14. Thesystem of claim 13, wherein the at least one computer processor isfurther configured for extracting lung measurements from thesegmentations.
 15. The system of claim 13, wherein the at least onecomputer processor is further configured for generating a disease mapfrom a visual comparison of the first series of images and the deformedimage.
 16. The system of claim 9, wherein the at least one computerprocessor is further configured for pre-processing the image data by oneor more of resizing to a standard resolution, scaling voxelattenuations, and affine-registering the second series to the firstseries.